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Abstract. Results of realistic simulations of solar surface convection on the 
scale of supergranules (96 Mm wide by 20 Mm deep) are presented. The simu- 
lations cover only 10% of the geometric depth of the solar convection zone, but 
half its pressure scale heights. They include the hydrogen, first and most of the 
second helium ionization zones. The horizontal velocity spectrum is a power 
law and the horizontal size of the dominant convective cells increases with in- 
creasing depth. Convection is driven by buoyancy work which is largest close to 
the surface, but significant over the entire domain. Close to the surface buoy- 
ancy driving is balanced by the divergence of the kinetic energy flux, but deeper 
down it is balanced by dissipation. The damping length of the turbulent kinetic 
energy is 4 pressure scale heights. The mass mixing length is 1.8 scale heights. 
Two thirds of the area is upflowing fluid except very close to the surface. The 
internal (ionization) energy flux is the largest contributor to the convective flux 
for temperatures less than 40,000 K and the thermal energy flux is the largest 
contributor at higher temperatures. This data set is useful for validating local 
helioseismic inversion methods. Sixteen hours of data are available as four hour 
averages, with two hour cadence, at steinr.msu.edu/~bob/96averages, as idl 
save files. The variables stored are the density, temperature, sound speed, and 
three velocity components. In addition, the three velocity components at 200 
km above mean continuum optical depth unity are available at 30 sec. cadence. 



1. The Simulation 

Solar surface convection on supergranular scales (96 Mm wide by 20 Mm deep) 
was simulated by solving the conservation equations for mass, momentum and 
internal energy. Spatial derivatives were ev aluated using sixth order finite dif- 



ferences (jNagaraian. Lele &: Ferziger 112003) and the time advance was by a low 



memory, third order Runge-Kutta scheme ( Kennedy. Carpenter &: Lewis 19991 ) . 
The calculations were performed on a grid of 1000^ x 500 giving a resolution of 
96 km horizontally and 12-70 km vertically, f-plane rotation is included corre- 
sponding to a latitude of 30. 

A tabular equation of state which includes local thermodynamic equilib- 
rium (LTE) ionization of the abundant elements as well as hydrogen molecule 
formation, was used to obtain the pressure and temperature as a function of 
log density and internal energy per unit mass. The radiative heating/cooling 
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Figure 1. Pressure, temperature, density, entropy (left), Gammal (solid), 
H ionization fraction (dashed) , and He ionization fractions (dash dot and dast 
dot dot dot) (right). The computational domain covers 10 % of the geometric 
depth of the solar convection zone and half its scale heights. There are 5 orders 
of magnitude in pressure, four in density and one in temperature within the 
convective region included. 



was obtained by solving the radiation transfer equation in both continua and 
lines using the Feautrier method, assuming Local Thermodynamic Equilibrium 
(LTE). The number of wavelengths for which the transfer equation is solved 
is drastically reduced by using a multi-group method whereby the opacity at 
each wavelength is placed into one of four bins according to its magnitude and 
the s ource function is binned the same way (jNordlund lll983 : IStein fc Nordlund ' 
20031 ). 



Horizontal boundary conditions are periodic, while top and bottom bound- 
ary conditions are open. The code is stabilized by diffusion in the momentum 
and energy equations, using a variable diffusivity that depends on the sound 
speed, the fluid velocity and the compression. 



2. The Atmosphere 

The computational domain covers five orders of magnitude in pressure in the 
convective region, which is half the scale heights in the entire solar convection 
zone, even though the domain is only 10% of the zones geometric depth (Fig. 
[1]). This includes the hydrogen ionization layer and the first and most of the 
second helium ionization layers. A shear layer develops with a shear of 100 m/s 
across the 20 Mm depth of the domain. 

Convective motions are turbulent. Descending vortex rings have typical 
mushroom shape and are connected back to surface with multiple twisted vortex 
tubes (Fig. [2]) . The horizontal scale of the convection increases from granule size 
at the surface to supergranule size near the bottom of the computational domain 
(Fig. E]). The horizontal velocity spectrum at the surface is a power law, with 
the magnitude of the velocity decreasing linearly with increasing size. There is 
a peak at granule scales and then a fall off at still smaller scales. About 2/3 of 
the area is upflows and 1/3 downflows nearly independent of depth. 
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Figure 2. Finite time (11.75 hours) Lyapunov exponent field ov a subdo- 
main 21 Mm wide by 19 Mm high, by 0.5 Mm thick, from a 48 wide x 20 Mm 
deep simulation. This corresponds closely to the magnitude of the vorticity. 
Figure courtesy Bryan Green, AMTI/NASA. 




Figure 3. Vertical velocity on horizontal slices at the surface and 2, 4, 8, 
12, and 16 Mm below the surface. Red and yellow are downflows. Blue and 
green are upflows. The dominant horizontal scale of the convection increases 
monotonically with increasing depth. 
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Figure 4. Contributions to the kinetic energy equation. Most of the buoy- 
ancy driving occurs just below the surface where the entropy fluctuations 
are large, but remains significant and nearly constant at larger depths. The 
buoyant work is balanced by the divergence of the kinetic energy flux near the 
surface and by dissipation in the interior. The change in behavior near the 
bottom boundary is due to boundary condition effects and is not physical. 

3. Driving and Damping 

Convection on the Sun is driven primarily by radiative cooling in a thin boundary 
layer at the bottom of the photosphere where radiation begins to escape to space. 
Radiative cooling produces low entropy, over dense fluid that is pulled down by 
gravity. There is also driving from the bottom of the convection zone where 
absorption of the radiative flux heats the fluid producing high entropy, under 
dense fluid. However, the entropy fluctuations near the surface are much much 
larger than those near the bottom of the convection zone, so most of the buoyant 
work occurs near the surface (Fig. . The buoyancy driving is balanced by the 
divergence of the kinetic energy flux close to the surface and by dissipation at 
small scales through the remainder of the convection zone. The damping length 
is 4 pressure scale heights (Fig. H]). 

4. Energy Fluxes 

Most of the energy in the cooler part of the convection zone close to the surface 
where hydrogen is partially ionized is carried by the internal (ionization) energy 
flux. At temperatures above 40,000 K, this decreases and most of the energy is 
transported by the thermal energy flux (Fig. [5]). 

5. Mixing Length 

The mass mixing length is the inverse of the logarithmic derivative of the unidi- 
rectional (either up or down) mass flux, I = {dF^s^ssidown) / dr / F^.^sidown))~^ . 
Figure [5] shows that over most of the computational domain except right near 
the surface the mass mixing length is 1.8 pressure scale heights. This is an 
easy, unambiguous, calculation to make for other stars as well to determine the 
appropriate mixing length. 
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Figure 5. Energy fliaxes (left) and mass mixing length (right). Note nega- 
tive fluxes are upward energy transport. The kinetic energy flux is downward 
and approximately half the magnitude of the total other, upward, energy 
fluxes. Energy is transported primarily as the internal (ionization) energy 
below 40,000 K and primarily as the thermal energy at higher temperatures. 
Upflowing fluid turns over and is entrained in the downflows in order to con- 
serve mass with a length scale or mixing length of 1.8 pressure scale heights. 



6. Application to Local Helioseismology 

Sixteen hours of data are available as four hour averages, with two hour cadence, 
at steinr.msu.edu/~bob/96averages, as idl save files. The variables stored are 
the density, temperature, sound speed, three velocity components. In addition, 
the three velocity components at 200 km above mean continuum optical depth 
unity are available at 30 sec. cadence. 



7. Summary 

A realistic simulation of supergranulation scale convection, covering half the 
scale heights in the solar convection zone can be used to evaluate local helioseis- 
mic inversion proceedures. It has also allowed us to determine the mass mixing 
length (1.8 Hp), to show that most of the energy is transported as ionization 
energy below 40,000 K and as thermal energy above that temperature, and that 
the buoyancy driving is largest close to the surface, but significant throughout 
the simulated domain. In the interior, this driving is balanced primarily by dis- 
sipation at small scales and the damping length is 4 Hp. Four hour averages 
with two hour candence of the density, temperature, sound speed, and three ve- 
locity components are available as idl save files. In addition, the three velocity 
components at 200 km above mean continuum optical depth unity are available 
at 30 sec. cadence. 
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